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^-j| The emergence of glassy behavior of the random field Ising model (RFIM) is investigated using 

an extended mean-field theory approach. Using this formulation, systematic corrections to the 
standard Bragg- Williams theory can be incorporated, leading to the appearance of a glassy phase, 
in agreement with the results of the self-consistent screening theory of Mezard and Young. Our 
approach makes it also possible to obtain information about the low temperature behavior of this 
glassy phase. We present results showing that within our mean-field framework, the hysteresis and 
avalanche behavior of the RFIM is characterized by power-law distributions, in close analogy with 

C"| ' recent results obtained for mean-field spin glass models. 
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I. INTRODUCTION 

The random-field Ising model (RFIM)lil is one of the simplest models used to describe the frustration introduced 
by disorder in interacting many-body systems. Despite its simplicity, its behavior proved to be the source of much 
controversy, primarily due to the lack of reliable theoretical methods that one can use for such systems. Still, efforts 
to elucidate the basic properties of the RFIM continue to attract considerable attention, primarily because of its 
direct relevance to a numbep,of important physical problems. These include not only the behavior of diluted magnets 
in external magnetic fieldsp but also several aspects of electronic transport in disordered insulatorsB and systems 
near the metal-insulator transition.ofl In addition, the non-equilibrium behavior of the RFIM has been used to model 
the physics of hysteresis and avalanche behavior and the origin of self-organized criticalityQ. Finally, non-disordered 
models with infinitesimal random fields have been studiedQ in order to investigate self-generated glassy behaviom 
observed in systems such as supercooled liquidsS, or even under-doped cuprates.Ea 

The simplest effect of turning on a weak random field is the resulting depression of the critical temperature for 
uniform ordering, while for sufficiently strong randomness the ordered phase completely disappears. This behavior is 
apparent even in the simplest Bragg- Williams (BW) mean-field theory (MFT) , but understanding the disorder-induced 
modification of the relevant critical beha"\por proved much more difficult. Very early on, perturbative renormalization 
group (RG) results of Parisi and Sourlasllil (PS) suggested the existence of a "dimensional reduction" by which the 
random problem belongs to the same universality class as a clean one in two dimensions less. 

Unfortunately, the beautiful result of PS was found, to be in conflict not only with the heuristic argument of Imry 
and Majl3 but also with-the exact results of Imbrie,Ej and thus deemed incorrect. More recently, the origin of these 
discrepancies was tracedi-3 to the implicit assumption of PS that outside the ordered phase, a single thermodynamically 
stable state exist. This assumption may not be warranted if another phase transition, presumably of glassy character, 
would precede any uniform ordering, which would lead to the breakdown of dimensional reduction. Of course, such 
a glassy phase is not found in the naive BWMFT, so that more sophisticated theoretical schemes have to be used 
in order to identify the corresponding instability of the high temperature paramagnetic phase. Such a theooi was 
first formulated by Mezard and Young,Ej who utilized the self-consistent screening (SCS) approach of Bray,E3 and 
identified the glassy phase by carrying out a replica-symmetry breaking stability analysis. Similar results were obtained 
by numerically solving the mean-field equations for a fixed realization of disorder by Lancaster et al.,li3 confirming 
the existence of the glassy phase. Finally, Mezard and MonassonJiH and De Dominicis et al.fLS presented arguments 
that the glass phase should persist even at weak disorder, and everywhere precede the uniform ordering, in agreement 
with the breakdown of the dimensional reduction. 

While these approaches provided important new information on the RFIM, several aspects remained unsatisfactory. 
The SCS scheme, while being able to identify the glass phase, proved of considerable complexity, making it difficult to 
obtain simple analytical results. In addition, the physical content of this formulation does not appear very transparent, 
making it difficult to establish what crucial ingredients are necessary for a theory in order to be able to identify and 
describe such a glass phase. Finally, the low temperature properties of the glass phase also seem very difficult to 
establish using this approach, which presents a severe linaitation in applying the RFIM to the problem of the "Coulomb 
glass"El and the related theory of Efros and Shklovskii.B 

The main goal of the present paper is to identify the simplest possible approach that is capable of providing a 
description of the glassy phase. We will show that to do this, two requirements have to be satisfied: (1) the theory 
has to identify the correct order parameters, and (2) it has to incorporate spatial fluctuations beyond the naive BW 
theory. Both of these are automatically satisfied by the SCS theory, but we will show how the same goals can be 
achieved in a simpler and more transparent fashion, by examining systematic corrections to the BW theory. To 
do thiSj-we follow an approach used by Plefkatj to derive the Thouless-Anderson-Palmer (TAP) equations for spin 
glasses.E3 A similar formulation was subsequently used by Georges et al.f3 to obtain systematic corrections to the 
MFT of the Sherrington-Kirkpatrick (SK) model. This approach fixes the desired order parameters by introducing 
appropriate source fields, and then evaluates the corresponding Gibbs free energy by an expansion in powers of the 
interaction J. A brief description of some of our results in the electronic context has already been presented in Ref. 
[5], but the present paper provides many more details and a number of new results. 

The rest of the paper is organized as follows. In section II, we begin our discussion by examining the RFIM on a 
Bethe lattice, where a particularly simple derivation of the BW theory and its leading corrections can be obtained, 
resulting in the emergence of the glass phase. A more general strategy based on the Plefka approach is presented 
in section III, showing that these leading corrections take the same form on arbitrary lattices. The low temperature 
structure of the glass phase is discussed in section IV, where the hysteresis and avalanche behavior in our extended 
MFTJs examined, showing the emergence of self-organized criticality similar to that recently discovered by Pazmandi 
et al.E3 for the SK model. The role of higher order corrections to our theory and its relation to the 1/z expansion 
approach is discussed in section V, where we show that the leading nontrivial corrections to the BWMFT represent 



dominant contributions in the joint limit of large coordination and large random fields. In this section, we also 
comment on the relevance of our approach to controversial question of self-generated glassy behavior in systems 
without disorder, which we discuss by examining the limit of weak random fields. Our conclusions are summarized 
in section VI, where we also outline some interesting directions for future work. 

II. THE BETHE LATTICE 

The simplest theoretical formulation of the random field Ising (RFIM) is obtained by considering the large coor- 
dination limit, where the standard Bragg- Williams mean-field theory (BWMFT) becomes exact. However, one has 
to go beyond this limit in order to obtain nontrivial results, and we address this question by examining systematic 
corrections to the BWMFT. The large coordination limit and the leading corrections are particularly easily formulated 
in the special case of the Bethe lattice, where a simple recursive procedure can be used. This approach also gives 
some insight in the mechanism for the emergence of the glass phase, so we begin our discussion by concentrating on 
the Bethe lattice in the limit of large coordination. 

The Hamiltonian of the random-field Ising model is given by 

Hint = — / JijSiSj — y ^ hiSi. (1) 

(v) 

Here, Si — ±1, and Jij — — are uniform ferromagnetic interactions between nearest neighbor sites, rescaled with 
the coordination number z in order to obtain a finite result in the z — > oo limit. The random magnetic fields hi are 
assumed to be Gaussian distributed, with zero mean and a variance (/if) = H RF . 

Using standard replica methodsEj, we can formally average over disorder, and the resulting partition function takes 
the form [a. = 1, ...,n; n — > 0) 

Z n = Tr exp[^ J2 £ S ? S ? + \(P H w? ££ $?)% ( 2 ) 

a fyj) i a 

We proceed by taking advantage of the tree-like structure of the Bethe lattice, by formally summing over all the 
degrees of freedom in one branch. The resulting functional &(Sq) is a function only of the variable Sff at the branch 
origin, and can be easily seen to obey the following self-consistent equation 

*(£?) = Tr s? (eM^-J2 S o S ? + ^W^E^f) 2 ]**- 1 ^)). (3) 



A. Bragg- Williams Theory 

In order to examine the large z limit, it is convenient to define a single-site effective action by 

L[S a ] = -±(f3H RF ) 2 (J2 S?) 2 - lnOF- 1 ^)). (4) 

a 

In the z — > co limit, this expression simplifies, since the interaction has been scaled by 1/z, and the functional <&(Sq) 
can be obtained by expanding the self-consistency condition, Eq. (3) in powers of the interaction J. To leading order, 
only the terms linear in J survive, and we find 

L^[S a ] = -l3J^S a m a - (PH RF ) 2 J2 sasfi - ( 5 ) 

a a <f3 

Here, the index (1) indicates that only terms linear in J are retained. As we expect for a Bragg- Williams theory, this 
local effective action corresponds to a single spin which, in addition to the local random field, also experiences the 
presence of a "molecular" field Jm a due to interaction with the neighbors. 

We emphasize that this procedure automatically defines the order parameters entering the local effective action. 
In the z — > oo limit only the magnetization m a appears, which from Eq. (3) satisfies the following self-consistency 
condition 



m a = (S a ) L[s] . 



(6) 



Since the interaction term does not mix the replicas, they trivially decouple, and the self-consistency condition for 
the magnetization m a = m takes the form 



m= Dxt&nh(((3H]iFX + j3Jm), 



(7) 



where Dx = dx exp(— x 2 /2)/\/2tt. This equation is a straightforward generalization of the well-known Bragg- Williams 
condition to include the effect of random fields. The critical temperature where the magnetization vanishes is easily 
computed, and is found to vanish at a critical strength H c RF j J = \/2/n of the random fields, as shown in Figure 1(a). 
Outside this ferromagnetic phase the local effective action of Eq. (5) reduces to that of noninteracting spins in the 
presence random fields, and no further phase transition can be found. 
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FIG. 1. Phase diagram for coordination numbers z — oo (a) and z — 4 (b). Ferromagnetic (FM), paramagnetic (PM) and 
glass phases are found. Note that the glass phase does not exist for z — oo. 



The reason for this limitation of the BWMFT is very simple, and can be appreciated by considering of the (thermally 
averaged) Weiss field h^ =< ^ JijSj >t produced by the neighboring spins on a given site. In a uniform system, 
1%Y is the same on every site, but in presence of randomness, it may display appreciable spatial fluctuations. This is 
especially important outside the any uniformly ordered phase, where the spatial average h =< lif > site, vanishes, but 
the local value hY of the Weiss field may remain finite, reflecting the local breaking of the up-down symmetry. This 

behavior is encountered in spin glasses, where fif = 0, but (hY) 2 becomes finite below a temperature corresponding 
to the glassy freezing. In the simple BWMFT, only the first moment hY ~ Ti is retained, thus its inability to describe 
any glassy ordering. 



B. Leading corrections and glassy freezing 



In order to search for the existence of nontrivial behavior outside the uniformly ordered phase, we have to go beyond 
the BWMFT, i. e. the z = oo limit. The leading correction is obtained by iterating the self-consistency condition Eq. 
(3) to second order in the interaction J, and in the expression for the local effective action a new term, quadratic in 
J appears 
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(8) 



a</3 



The order parameter m a is still given by Eq. (6) (although the average is now computed with respect to the extended 
action of Eq. (8)), but the new order parameter q a @ appears, which is self-consistently determined by 



n a/3 



(S a S ) 



L[S]- 



(9) 



This quantity is nothing but the familiar Edwards-Anderson order parametral 2 ^. Note that in our Bethe lattice 
procedure, it appears automatically, as a result of an expansion in powers of the interaction J to the lowest nontrivial 
order beyond the BWMFT. Its presence reflects the fact that in finite-coordinated lattices, local Weiss fields hi are 
random numbers with a finite dispersion. More precisely, they are described by the distribution of local fields P(hi), 
the explicit form of which is determined by the replica matrix q a/3 . 

As in standard spin-glass theoryc3, the solution of these equations assumes the simplest form in the high temperature 
phase, where replica symmetry is valid, such that m a = m and q a @ = q. In this case, our self-consistency conditions 
take the form 



m = Dx tanh 
q = Dx tanh 2 



(0H RF ) 2 + ^-(q-m 2 )\ x + fiJm 
(PH RF f + &-{q-m 2 )\ x + pjm 



(10) 



This is sufficient to determine the ferromagnetic (FM) phase boundary which, as before, is determined by setting 
to = 0. By numerically solving these equations, we find that the FM phase is only slightly reduced due to the 
fluctuation corrections, as shown in Figure 1(b). 

Identifying the glassy freezing is more difficult. For standard spin- glass models, the glassy freezing coincides with 
the breaking of the up-down symmetry, so that the glass transition temperature can be identified even within replica 
symmetric theory, as the point where the Edwards- Anderson order parameter q — < Si > 2 assumes a nonzero value. 
In our case, the random magnetic field plays a role of a source conjugate to the order parameter, locally breaking 
the up-down symmetry. The situation is similar as in spin glass models in a uniform external field 2 !], where the 
replica symmetric order parameter q remains nougjero for any temperature, and thus cannot be used to identify glassy 
freezing. Instead, we follow Mezard and YoungJiJ, and look for an instability to replica symmetry breaking (RSB) 
within the paramagnetic phase. To do this, we can set to = 0, aad note that the remaining equation for g Q/3 is in fact 
identical to that describing the Sherrington-Kirkpatrick mode£3 in presence of random magnetic fields. This model 
is also described with the Hamiltonian of Eq. (1), but this time with J^-'s being Gaussian random variables with zero 

mean and variance («/£■) = ■jfjj, where J is the interaction of the original lattice model. 

The advantage of mapping our equations to those of an infinite range model will be also used in sections IV and 
V, where we examine the low temperature structure of the glass phase. In addition, the .calculation for obtaining the 
replica symmetry breaking (RSB) boundary can be carried out using standard methodscJCJ, following the approach 
of De Almeida and Thoulesscd (AT). A more general strategy for performing the RSB stability analysis, which is 
applicable for arbitrary lattices, and to higher order fluctuation corrections will be presented and discussed in section 
V. Here we just quote the result valid to the leading nontrivial order, which takes the form 



(/3J); I Dxcosh- 4 



(/3H 
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{(Mfq 



1/2 



(11) 



As expected, in the large coordination (z — > oo 
reduce to standard BWMFT predictions. For finite z 



limit, the glass transition temperature vanishes, and our results 
the SG phase emerges, in agreement with the results of Mezard 



and Younglil. The above equation can easily be solved numerically, and the results are shown in Fig. 1(b). In contrast 
to the predictions of the SCS approach, our SG phase emerges only for sufficiently strong random fields. This is a 
simple result of the fact that in simple mean-field treatments such as ours, different phases do not directly affect each 
other, since there is no "precursor" of the ordering in the disordered phase. 



C. Limit of large random fields 

It is interesting to note that the situation is especially simple in the limit of large random fields. In the extreme 
case Hrf — > oo, the interactions can be ignored and all the spins tend to align with their local random fields. In 
this case there is only one thermodynamic state, so that the emergence of a multitude of metastable states associated 
with the RSB instability is clearly suppressed. We conclude that the glass transition temperature must be depressed 
to zero as Hrf —> oo, just as in the case of spin glasses in uniform external fields first discussed by de Almeida and 
Thouless. However, the asymptotic form of Tg(Hrf) proves to be different in our case, which may be significant for 
several experimental systems. 

In the limit of large random fields, our results simplify considerably, since the RSB boundary resides at very low 
temperatures. The leading order behavior can be obtained by replacing the q(T) in Eq. (11) by its zero temperature 
limit q(0) = 1, and we find 



T G 



4J 2 



ZzH 



RF 



1/H RF . 



(12) 



It is interesting that our glass transition temperature ftps decreases very slowly with the random field strength, 
in contrast to the form of the de Almeida-Thouless lineEaiHj(infinite ranged spin glasses in a uniform field) where 
Tg ~ exp(— H 2 /2J 2 ). This fact could be particularly significant for strongly disordered electronic systemsQ'U, where 
it would suggest the possibility to observe the glassy behavior of electrons at finite temperatures. The physical reason 
for this behavior in our case is not obvious, but we will see that it reflects some very subtle features of the low 
temperature glass phase, which will be discussed in section IV. .-, .— . 

From the physical point of view, an instability to RSB, such as we find, is knownErES to describe the emergence 
of an extensive number of metastable states, and the associated slowing down in the relaxational dynamics of the 
spins. Experimentally, this results in the onset of the history dependence of cooling, and the related bifurcation of 
field-cooled (FC) and zero-field cooled (ZFC) spin susceptibilities. 
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FIG. 2. Glass transition temperature as a function of the random field strength. Open circles are the experimental data for 
Feo.31.Zno.69Fo. 2 (from Ref. [2]), and the line is the prediction of Eq. (11). The inset shows how Tq ~ 1/ Hrf at large fields. 



Such experiments have been carried out on diluted antiferromagnets in uniform external fields, which have long been 
believed to be realizations of the RFIM. Here, the effective strength,of the random field can be varied by modifying the 
magnitude of the external magnetic field. In one such experiment^ the field dependence of the "irreversibility line" 



has been determined, defined as the temperature where the FC and ZFC susceptibilities start to differ. Interestingly, 
these experiments show a rather slow decrease of this glass transition temperature, in agreement with our predictions. 
We have digitized the data from Ref. [2], and compared them to our predictions, as shown in Fig. 2, where an 
apparent confirmation of our Tq ~ 1/ ' Hrf law can be seen (see inset). While this agreement of our theory and 
experiment is encouraging, it should not be taken too seriously, given the uncertainties of the precise correspondence 
of the experimental system and the RFIM that we consider. More experiments on other related systems would be 
welcome to test our predictions in more detail. 

III. GENERAL LATTICES: LEGENDRE TRANSFORM APPROACH 

So far, we have seen how leading corrections to the BW theory can be obtained on the Bethe lattice, resulting in 
the emergence of the glass phase. This example was useful, because it automatically introduces the correct order 
parameters, thus allowing the emergence of the glass phase. However, we would like to formulate a more general 
approach, in order to demonstrate the generality of our conclusions, and also to be able to systematically examine 
higher order fluctuation corrections. In the case of the RFIM, the desired order parameters cannot be introduced as for 
standard spin-glass models, where one decouples the disorder-averaged interaction term using a Hubbard-Stratonovich 
transformation. In absence of random interactions, the "bare" disorder-averaged Hamiltonian has only terms linear 
in the interaction J, but higher order terms can be generated by fluctuations. In this case, one has to introduce the 
order parameters "by hand" , and then expand the free energy to the lowest nontrivLal order in the interaction J, in 
order to obtain glassy ordering. To do this, we follow a strategy introduced by Plefkalij and Georges et al.Eil, and use 
a Legendre transform approach, introducing external source fields £" in order to fix the Edwards- Anderson (EA) 

order parameters qf . Here and in the rest of the paper, we are mostly interested in the emergence of the glassy 
phase of the RFIM, so we concentrate on the non-magnetic (to = 0) solutions. The disorder-averaged Hclmholtz free 
energy takes the form 

-0F= lim A [-0F n ] , (13) 

n->0 an 

where 

- PF n = lnTrexp[^ £ £ S?S? + V*f£ £(£ Sff + £ £ S*S?S?], (14) 

a (-ij) i a z q</3 

The EA order parameters are given by 



■ -.■_ d(-0F) ,.,„.,., 



^f = ^£^ = (S?Sf), (15) 



and the corresponding Gibbs free energy is 



(3G = lim A [-f3G n ] ; -0G n = -[3F n - E E tfV- (16) 

i a</3 



A. Extended mean-field theory 

To obtain an extended mean-field theory, we can expand — (3G to a desired order in powers of the reduced interaction 
e = f3J, while fixing the value of the independent variables q" . In doing this, one has to keep in mind that source 
fields £" are by Eq. (15) implicit functions of the order parameters qf . The form of these functions is also a function 
of PJ, so in order to be consistent, one has to expand the fields £ a @ to a given order in (3 J as well, while fixing qf-s. 
Defining —[3G n /N = g, we can write 

.9(e) = 9o+gtnt(e), (17) 

00 e k 
9int{e) =J2j\9k- (18) 

fe=i 



In this expression, the coefficients gk are functions of q a @ evaluated at e = 0. Explicitly, we find 



»m=/i]4eecv. as) 



where 



TV 
and the reference effective action is defined by 



9u^ : ,..,., N 

i a</3 



/ [£] = 4lnTrexp{-L }, (20) 



Lo = -Y,H {tf + (PH RF f) SfSf . (21) 

i a<j3 



The interaction terms take the form 



.91 = ^{4>)o (22) 

92 = ^(<K<l>-{<l>)o + X-(x)o))o, (23) 



where we have defined 

^=\Y.Y. S ? S t' ( 24 ) 

X = EE^- (25) 

i a<P 

In these expressions, the averages are taken with respect to the reference effective action L 0l which is a function of 
the external fields £" . Note that since the coefficients gk are evaluated at e = 0, after the cumulants are evaluated, 
the external fields should be set equal to £™ (e = 0) = £"^, which are implicit functions of the independent variables 
q , as determined by the following constraint condition 



r) 

a af3 = f \f] 



_ I not cfi 



(S?S?) . (26) 



tf"=tf" 



At this stage, we have taken the order parameter q a @ to be uniform, since the system is translationally invariant after 
averaging over disorder. 

The equation of state is obtained from the saddle-point condition, giving 

= ffS = E^/o - Q-ys]d a0 Q S - C P + d a p 9int, (27) 

7 5 

where we have used the notation d a p = d/dq a °, and d a p = d/d£ al3 . Using the constraint condition Eq. (26), we can 
eliminate £"^ and write 

q a0 = d a0 f o [d g int }. (28) 

The described expansion can be carried out to any desired order in e, generating fluctuation corrections to the Gibbs 
potential. It is worth noting that, since this expansion is a series in powers of the interaction Jy, it generates 
corrections due to short-ranged fluctuations, in contrast to the usual loop expansion that is dominated by the long- 
wavelength fluctuations. Accordingly, we expect that the prediction for the critical behavior of the order parameter 
and response functions to maintain the mean-field character if the expansion is truncated at any finite order. 



B. Leading fluctuation corrections: the J 2 theory 

The leading fluctuation corrections to the BW theory are obtain by retaining terms up to second order in e, which 
is the lowest order to which we can identify a glass phase. Using the up-down symmetry of the Hamiltonian, we find 
(Appendix A) that {<fi} = and \{ £ = 0) = 0, so that the above expressions simplify giving 

.91 = (29) 

.92 = ^(0 2 >o = ^[n + 2^ (7 Q V /3 ]. (30) 

a<{3 

This expression is valid to any order in e, but using in the leading theory (terms to 0(e 2 )), we find 

e 2 
d a pgvnt = — g a/3 , (31) 

z 

giving 

q aP = d af} f [- z (f]J) 2 q]. (32) 

This equation of state is precisely the same as the m = limit of the expression (Eq. (9)) that we have obtained on 
the Bethe lattice. We have thus demonstrated that to leading nontrivial order beyond the BW theory, we find the 
same glassy phase for arbitrary lattices. To this order, which corresponds to introducing the Onsager reaction field 
correction, the results are independent of the lattice geometry, as also found in other problemsEJEi 

IV. THE T = GLASS PHASE: SELF-ORGANIZED CRITICALITY 

The nature of the glass phase is in general much more complicated than the usual ordered phases, and is characterized 
by a large density of low energy excitations and low-lying metastable states. Within the class of mean-field models 
for spin glasses, the structure of the ordered phase has been investigated most extensively in the vicinity of the glass 
transition, where the Parisi theory is analytically tractable. Much less is known about the low temperature behavior, 
although some of the most interesting phenomena are most pronounced there, including the hysteresis and avalanche 
behavior. 

A. Hysteresis and avalanches behavior in presence of random fields 

To investigate the glass phase, we concentrate on the T = behavior of our model. To do this, we limit our attention 
to the J 2 theory described above, where the RFIM maps to the SK spin-glass model in presence of additional random 
magnetic fields. This mapping allows us to use different methods to investigate the T = behavior, which would 
be difficult to address using the Parisi theory. In particular, we investigate the hysteresis and avalanche behavior 
within— our extended mean-field theory of the RFIM, by following recent work of Pazmandi, Zarand and Zimanyi 
(PZZ)E 2 ! In their original calculation, PZZ have examined the standard SK model, and the only modification that we 
introduce is the additional presence of random magnetic fields. The procedure that we follow is as follows. One first 
introduces a large external magnetic field (in addition to fixed random fields), so that all the N spins are aligned with 
it. The external field is then slowly reduced and the stability of the spin configuration is examined with respect to any 
single spin flips. As soon as the system becomes unstable, the spin configuration is allowed to relax to a local energy 
minimum, causing an "avalanche" . The procedure is then repeated, resulting in the system following a hysteresis 
loop. If the external field is sweeped to large negative values (such that all the spins align in the negative direction), 
and then the procedure reversed, then the state of the system follows a "major" (external) hysteresis loop. If the field 
is instead reversed before the major loop is completed, then the system embarks on one of the minor hysteresis loops. 
This procedure is illustrated is Fig. 3 where a typical "hysteresis spiral" is presented. 




FIG. 3. A typical hysteresis spiral. 



B. Distribution of local fields 



To characterize spin- glass state, we examine the probability distribution P{Hi) of local magnetic fields Hi = 
hi + y\ JijSj acting on a given spin S%. In the high temperature phase, P(Hi) has a simple Gaussian distribution. 
To see this, we recall that above the glass transition, the replica symmetry remains valid, so that the presence of the 
interaction term in the local effective action (see Eqs. (9-10); we use m = here) simply renormalizes the effective 
distribution of random fields, such that H^J, = \/(Hrf) 2 + J 2 q/z. The distribution P(Hi) therefore remains a 
simple Gaussian of width given by H^J, . From Parisi theorycEl, we expect this distribution to acquire a nontrivial 
form upon replica symmetry breaking, but we would like obtain its specific form in the spin-glass state. 

To calculate this quantity on the hysteresis loop, we,follow a simulation procedure identical to that used for the 
SK model by Pazmandi, Zarand, and Zimanyi (PZZ)E3. In this procedure, one considers a finite size system of N 
spins, with a given realization of random interactions Jy and random fields hi, and let the system explore the meta- 
stable states sampled on the hysteresys loop. The values of the local fields Hi are then computed for every spin, 
and the procedure is repeated for many realizations of disorder in order to generate a large ensemble from which the 
desired distribution can be computed. To implement this procedure, we have carried large-scale simulations using 
systems with up to N = 3200 spins, and obtaining ensembles of M = 500, 000 data points from which the distribution 
histograms were obtained. The resulting distributions for several values of the random field strength are shown in 
Fig. 4. We find that the distribution is characterized by the emergence of a universal pseudogap of the form 



P(Hi)^CHf, {Hi«J) 



(33) 



where ct-j= 1 and C = z/J 2 is independent of the random field strength. This universality is consistent with the findings 
of PZZl3, who have shown that it reflects the self-organized criticality which characterizes the T — behavior of the 
mean-field glassy systems. Our results confirm that the conclusions of PZZ remain valid in presence of random fields, 
thus representing a very robust property of glassy phases, at least within the confines of the considered mean-field 
descriptions. 
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FIG. 4. Distribution of local fields for N=3200, as a function of random field strength^ for Hrf^Jz/ J 
universal form of the pseudo-gap. The finite value of P(0) ~ 1/y/N is a finite size effec1c3. 
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By following these procedures, we have carefully verified that all the findings that PZZ have established^ for the SK 
model also hold in presence of random fields, thus also apply to the RFIM within the present mean-field formulation. 
In particular, we have confirmed that the avalanche sizes are distributed on all scales, and are characterized by a 
power-law distribution, characteristic of self-organized criticality. It is most remarkable that this critical nature is not 
confined to the ground state, but persists for all the metastable states within the hysteresis loop. It is particularly 
interesting, as we have explicitly verified by simulation, that not all local minima of the energy have this property. 
Instead, the critical states form a subset of metastable states that can be reached by the described hysteresis procedure. 
In this way, the ground state seems not to have any special features, but rather to share the same properties with 
an extensive number of critical metastable states. This notion offers a natural origin for the criticality that is found, 
since all the states along the hysteresis loop can be considered to be on the brink of an avalanche, and are therefore 
inherently unstable to weak perturbation, making the criticality possible. In fact, it is precisely the requirement 
for marginal stability of the hysteresis states that was used by PZZE3 to derive the universal form of the local field 
distribution. This argument examines the modification of the local fields Hi upon flipping a given set of nfu p other 
spins. This is given by 



6Hi = H' - Hi 



2 £ 

j flipped 



Jij bj . 



(34) 



One then computes the probability that the local field is reversed, so that instabilities are created, triggering 
avalanches. Note, however, that the variations SHi is independent of the value of the (external) random fields hi 
present in our case. As a result, the rest of the argument goes as in PZZE3, giving the above marginality condition. 

From the historical perspective, evidence of marginal stability of the spin glass_phase in mean-field models has long 
been appreciated based on Parisi and Thouless-Anderson-P-almer (TAP) theoryEI In addition, stability arguments 
requiring a > 1 have been presented by Palmer and PondEa, based on early ideas of TAPES. However, it was not 
clear why this bound has to be satisfied, or that even the prefactor C assumes a universal value. In this sense, the 
hysteresis and avalanche study of PZZ has provided an important conceptual advance, linking this universality with 
the self-organized criticality of the mean-field models. 

C. Self-organized criticality and the AT line 

We have seen that the glass transition temperature in our random- field case decreases very slowly, as Tq ~ 1/Hrf 
at large random fields. This is very different than in the case of the familiar AT lineEj of the SK model in a uniform 
field, where the RSB temperature decreases exponentially with the uniform field strength. In the following we present 



11 



a simple heuristic argument that explains the physical origin of both behaviors based on the self-organized criticality 
of the glassy ground states in these mean-field models. 

What emerges from the analysis of PZZ for the SK model, and which also applies in our extended mean-field 
formulation for the RFIM, is the phenomenon that at T = the distribution of local fields assumes a universal 
form at Hi — > 0, as given by Eq. (33). We have established this property for all the states within the hysteresis 
loop following the methods of PZZ. However, it is very likely that a similar universal distribution of local fields also 
characterizes the exact ground statc_of the system. In fact, this possibility has been proposed a long time ago both 
by TAPE3 and by Palmer and PoncO, consistent with the notion of the marginal stability of the spin glass state for 
infinite-range models. Similarly to those of PZZ, the arguments of Palmer and Pond can also be extended to include 
the addition of random magnetic fields, resulting in a stability bound that remains universal, i. e. independent of the 
random field strength. To establish more firmly this universal form of the pseudogap for the ground state, one would 
have to carry-out more elaborate numerical simulations. Efficient optimization methods needed for such computations 
are available^, but this requires extensive efforts which are beyond the scope of this paper. 

For our purposes, we will assume that the pseudogap retains a universal form for all T = critical states, including 
the hysteresis states as well as the ground state. If this is true, we can now make a simple estimate of the "condensation 
energy" gained by glassy ordering, which should scale with the gap size. Having in mind that the replica symmetric 
(high temperature) distribution P_Rs(Hi) is a Gaussian of width given by H R p , and the fact that the pseudo-gap has 
a universal, linear form, we conclude that the gap size should scale as 

E gap ~ Pbs(0) ~ 1/Hmr, (35) 

at Hrf — > oo. Here, we have used the fact that q(T = 0) = 1, so to leading order H R p sa Hrf for large random 
fields. Using this result, we can immediately estimate the glass transition temperature, which should also scale with 
this gap size, since E gap is the only energy scale characterizing the ground state. We thus conclude that T g ~ 1/Hrf, 
in agreement with our analytical calculations based on the RSB analysis. It is worth emphasizing here that the 
simple relationship between the gap energy and the random field strength directly follows from the universality of 
the pseudo-gap form. Has the gap had another functional form (e. g. a randomness-dependent exponent a), this 
relationship would be modified, resulting in a different dependence of T g (H RF ). In this way, the agreement between 
out analytical RSB results and the presented heuristic arguments provides additional evidence for the self-organized 
critical nature of the ground state in these models. 

Finally, it is interesting to examine how the presented heuristic argument applies to the usual SK model in a uniform 
field. In this case, the only modification is introduced by the fact that in the high temperature phase, the uniform 
field simply shifts the entire Gaussian distribution of local fields. As a result, Prs(H{ — 0, H) ~ exp{— H 2 /2 J 2 }j-and 
we find the AT temperature decreases exponentially with uniform field, in agreement with the analytical results 2 ^ of 
AT. 

D. Relevance for the Coulomb glass problem 

The nontrivial nature of mean-field glass models may be particularly significant for the Coulomb glass problem, 
relevant for disordered insulators. The Hamiltonian for the Coulomb glass is given by 

e 2 
H = 22 —( n i - ( n ))( n j ~ (%•)) -^2 Sim. (36) 

Here, n,; = 0, 1 are the occupation numbers, ry is the inter-site distance, and £i are the random site energies. 
The transformation Si — 2rii — 1 immediately maps this Hamiltonian to an antiferromagnetic RFIM with long- 
range interactions. Note that the presence of long-range antiferromagnetic interactions leads to considerably stronger 
frustration than that in the ferromagnetic RFIM considered in this paper, which is found only for sufficiently strong 
randomness. These differences may be crucial in low dimensional systems, possibly leading to low temperature glassy 
phase for the Coulomb system even if a similar behavior is suppressed for the ferromagnetic RFIM case. However, on 
the mean-field level, even the standard RFIM displays displays such glassy ordering, the character of which may be 
closely related to that of the Coulomb glass. 

The best established property of the Coulomb glass model is the existence of "Coulomb gap" , as predicted by Efros 
and Shklovskii (ES)I3. According to ES, the electronic system would be unstable, unless the single-particle density of 
states (which corresponds to our local field distribution) has a pseudo-gap of the form 

p{e) = C{d)e d -\ (37) 
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where C(d) is a universal constant in d dimensions. More precisely, ES have examined the stability of the system 
with respect to one-electron excitations, showing that the form of Eq. (37) represents an upper bound for the density 
of states. If this bound would be saturated, then the pseudo-gap would assume a universal form, hut, no convincing 
arguments have been presented why this should happen. However, large-scale numerical studiescJ have obtained 
results similar to the ES predictions, failing to produce any evidence of the hard gap. Still, the reason for saturating 
the ES bound has remained a mystery. 

Another aspect of the Coulomb glass that has not been properly clarified is the presence or absence of a finite- 
temperature glass transitiono, and the nature of the glass phase. In this respect, the main difficulty was the absence 
of an appropriate order parameter that would allow identifying the transition. Since the random site energies in this 
model play a roie,of random fields, the situation is identical as in the usual RFIM, and the usual EA order-parameter 
cannot be usedc3. As we have seen, the transition at best can have a character of an AT line, which should be 
most easily identified in changes of the dynamics as the temperature is lowered. Some numerical eiddence that such 
a dynamical transition may exist has been obtained by using the "damage-spreading" algorithmst3, giving hints of 
crgodicity breaking below a certain temperature. In addition, several experiments^! have reported history-dependent 
transport and other glassy features in disordered insulators. 

Despite the successes of the ES theory, several basic questions remain unanswered. Most importantly, it is not clear 
whether the emergence of the universal ES gap is related to the possible low temperature glassy state of the model. 
In this respect, the scenario that we have presented for the RFIM provides an interesting possibility. It is conceivable 
that, as in our mean-field theory, a glassy phase exists below a well-defined transition temperature, which corresponds 
to the emergence of a large number of metastable states. If this phase were characterized by self-organized criticality, 
then the associated marginal stability would naturally explain the saturation of the ES bound, and the resulting 
universality of the Coulomb gap. An an interesting way to address these questions would consists of examining 
the T = hysteresis properties of this model, following the approaches of PZZ, but this direction will be pursued 
elsewhere. 

V. HIGHER ORDER CORRECTIONS AND THE 1/Z EXPANSION 

So far, we have examined the leading-order corrections to the BW theory, producing the glass phase. On general 
grounds, one does not expect that higher-order fluctuations corrections of a finite order in J would produce qualitative 
modifications in this mean-field context. Nevertheless, we shall examine the ncxt-to-leading terms, in order to assert 
the convergence properties of such an expansion. 

A. The J 4 theory 

In section III, we have already calculated the terms up the to order J 2 in the expansion of the Gibbs free energy. 
Terms up to order J can be computed using the same procedures, where we use the average up-down symmetry 
of the problem outside the ferromagnetic phase. In contrast to the J 2 theory, these additional terms depend on the 
specific form of the lattice. To be specific, we concentrate on the hypercubic lattice in d dimensions (so that z = 2d). 
After lengthy algebra, the resulting contributions up to order J take the form 

93 = 0, (38) 

.94 = ^[<0 4 >o + 3(^V>o]- (39) 

Here, \' = dx/de 1 and the brackets (• • -) c indicates that only connected diagrams, which give non-vanishing contri- 
butions, should be retained. This refers to diagrams obtained by representing the quantity <p (see Eq. (24)) by a 
bond, since it stems from the interaction terms connecting two nearest neighbors on a lattice. While each power of <£> 
involves a sum over all possible embeddings of such a bond, only diagrams consisting of a close loops of such bonds 
produce non-vanishing results, due to the average up-down symmetry of the problem. Four different classes of such 
diagrams are shown in Fig. 5, corresponding to the J contributions that we consider. 
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(a) (b) (c) (d) 

FIG. 5. Diagrams contributing to order J 4 . The diagrams (a), (b) and (c) correspond to the (4> 4 )o term, and (d) to the 
b 2 x'}o term. Here, the 4> bonds are represented by a full line, and the local x' term is represented by an open circle. 



The evaluation of these terms is straigthfbrward. Here, we only emphasize the following properties of these diagrams, 
that are valid even to higher order in J, as follows. (1) Diagrams (a), (b) and (d) contain vertices with four emerging 
bonds, that in the replica calculation give rise to expressions which include moments of the form r a pjs = (S"S^ S? Sf) , 
with a 7^ j3 7^ 7 7^ 5. Such moments cannot be simply expressed in terms of the order parameter q a @ . In contrast, 
the "loop" diagram (c) is expressed as a power series involving only powers of q a @ . (2) The lattice embedding factors, 
which specify the z-dependence, do not depend on the specific replica decorations that have to be carried out on 
each diagram, but are determined only by the topology of the diagram. (3) We find that diagrams (b) and (d) result 
in identical expressions, except for the lattice embedding prcfactors, which are different. As a result, cancellations 
occur, such that the sum of (b) and (d) produces a contribution to the free energy which is of order J 4 /z 3 , as is the 
contribution of (a). (4) The leading contribution, of order J 4 /z 2 , follows only from the loop diagram (c). (5) The 
diagrams (a), (b), and (d) consist of self-retracting paths, and as such are identical for both a hypercubic and the 
Bethe lattice with the same z. The "open loop" diagram (c) is specific to the hypercubic, but absent for the Bcthe 
lattice. 

In examining the J 4 (and higher order) corrections, we are primarily interested in identifying a limit where these 
terms are "small" , so that only the leading J 2 contributions may be retained. On a Bethe lattice, all the contributions 
are of order l/z 3 , thus can be ignored in the limit of large coordination. It is easy to see that the same coj&clusion 
applies to arbitrary lattices with purely random interactions such that (Jy) = 0, as found in spin glass modelscil, since 
in this case only self-retracting paths survive. In contrast, for the RFIM on general lattices, where the interactions are 
uniform, the "open loop" diagrams (c) survive, and provide the leading contribution, which is of the same order as in 
the J 2 term. We expect that similar conclusions apply to higher order contributions as pfll- This property, that the 
"open loop" diagrams provide a leading contribution in large dimensions, is well knownEil, and has been extensively 
used to study models of strong electronic correlation in large dimensions^. To get the leading contributions for large 
z, one would therefore have to sum up s . yic h open loop contributions to all orders in J. In absence of randomness such 
re-summations have been carried outEilS, but at least for the electronic models, the results were not qualitatively 
different from those obtained on the Bethe lattice, where only the second order terms may be retained. 

In our case such re-summations are more difficult to carry out due to the presence of randomness (replica indices), 
and will not be attempted here. Instead, we will show that if one examines the joint limit of large coordination and 
strong random fields, even the J 4 terms represent sub-leading contributions, and the simple J 2 theory suffices. To 
show this, we only consider the leading loop contributions (diagram (c)), which takes the form 
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B. 1/z corrections to the equation of state: RS solution 

A general expression for the equation of state, valid to arbitrary order is given by Eq. (28). To compute the relevant 
J corrections in the replica symmetric case, we calculate the variation of gi nt , which to this order reads 



d a p g%- 



I AS 



-q 



q(l-q) 2 , 



(41) 



where we have taken the n — > limit. 

We are especially interested in examining the form of the solution in the large Hrf limit, where the RSB transition 
occurs at low temperatures. We examine the relative magnitude of the terms appearing in gi n t- At first sight, the J 4 
terms seems to be dominant in the low T limit, since it is of order e 4 = (/3J) 4 . However, note that this term is also 
proportional to 8q 2 = (1 — q) 2 , which is small at low temperatures, since q — > 1 at T — ► 0. To see this, it suffices (to 
leading order) to compute q(T) at J = 0, giving 
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^=o = l_ 2 1 ... (42) 

y/ir pH R F 

This gives Sq ~ ((3H RF )-\ so that (/3J) 4 q(l - <?) 2 - {f3J) 2 {H RF / J)" 2 . The J 4 term is down by a factor {H RF /J)- 2 , 
thus contributing only to sub-leading order in the limit of large random fields. 

C. RSB stability analysis: a general approach 

Within the J 2 theory (sections II and III) , we have mapped the RFIM to an equivalent infinite-range model, making 
it possible to carry out the RSB analysis similarly as for the SK model. In the following, we present a general approach 
to the RSB stability analysis, which can be used-even if higher order terms are retained. 

To identify the RSB instability we follow ATEi, and examine the variations of the Gibbs free energy with respect 
to the deviations from replica symmetry q a/3 — q + Sq a ^ . From the definition of the Gibbs free energy, Eq. (19), the 
corresponding stability matrix is 

d a 0,y5 g — —dapC + d a p^& Qint- (43) 

Here, we have used the notation d a p^s = d 2 / dq a @ dq lS . To eliminate the quantity d a p(; lS , we take a variation of the 
constraint condition Eq. (26), and obtain 

= 6 a0 , lS - J2 Sr*,/* fo[£\d a p^\ (44) 

where we have used d a p^$ = d 2 /<9£ Q ^<9£ 7 * . To write these expressions in a more compact form, we introduce a matrix 
notation {g") a p. :1 s = d a p. n s g; {g'Lt)otP,iS = d a p. n s g m f, (£'W,7<5 - d <*p V b \ (/mtW, 7 <5 = d a/ 3 n s /<,[£], giving 

9" = -i' + 9int (45) 

= 1- /"•£'. (46) 

In this form, the matrix £' can be eliminated and we get 

g' = -{f"r l + gL- (47) 

To simplify the calculation further .-we note that in the high temperature phase all the eigenvalues of g" are positive, 
but the RSB instability is identifiedE3 when at least one of its eigenvalues vanishes, such that the determinant (at 
fixed n) vanishes. Using this property, an equivalent approach to the RSB stability analysis can be formulated by 
examining the stability of an auxiliary matrix 

g>( = f"-g". (48) 

This is true, since 

det(g'{)=det(f")det(g"), (49) 

and it can be verified by explicit calculation that the matrix /" remains nonsingular in the temperature range of 
interest. In other words, to identify the RSB instability, we need to compute the eigenvalues of the auxiliary stability 
matrix 

g" = f"-g'L-i- (so) 

We now present a general strategy for computing the eigenvalues of this matrix. Our first observation is that any 
replica matrix of the form considered, has at most three different matrix elements, as discussed by ATcij. In additions, 
for any such matrix there can be found three different (degenerate) eigenvectors. Most remarkably, the form, of these 
eigenvectors does not depend on the vale of the corresponding matrix elements. Thus, any replica matrix has the same 
eigenvectors, and only the corresponding eigenvalues depend on the value of the matrix elements. Since a product of 
two replica matrices is again a replica matrix, we conclude that the relevant eigenvectors in our case are identical to 
those computed by AT, and we only need to evaluate the relevant eigenvalues. Only one of the eigenvalues vanishes at 
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the RSB transition, call it A3, and the corresponding eigenvector X3. The eigenvalue A3 in our case can be computed 
by acting on X3 with the matrix g'{, and we find 

A3 - \{ ■ A 3 "* - 1, (51) 

where A3 and A™* are the respective eigenvalues of the matrices /" and g'/ nt . We emphasize that this strategy is not 
specific to our J 4 theory, but is valid to arbitrary order in the expansion. 

To obtain the desired J 4 corrections, we need to compute the corresponding corrections to the matrices /" and g" nt , 
as done in Appendix B. For general Hrf, we find that the J 4 corrections are of the same order as the J 2 ones, leading 
to a glass transition temperature Tq ~ l/\/z. Clearly, all the higher order terms coming from the "loop" diagrams 
are also of the same order, and would have to be included as well, in. order to collect all the leading contributions in 
a 1/z expansion, similarly as in other problems in large dimensionsE3. 

However, as shown in Appendix B the expressions simplify in the limit of large random fields, where all the J 4 
contributions are down by a factor (J/Hrf) 2 , and to leading order expressions are obtained by simply retaining 
only the J 2 terms. We conclude that the J 2 theory, which provides the leading nontrivial order in our extended 
mean-field theory, represents an asymptotically exact formulation in the joint limit of large coordination and large 
random fields. These conclusions have been obtained by examining the example of a simple hypercubic lattice with 
nearest neighbor interactions, in the limit of large coordination. Since the general structure of the diagrammatics and 
the relevant z-dependence is qualitatively the same for general lattices, we expect these results to hold for any model 
with short-range interactions. 

The situation is more complicated for models with long-range interactions, such as the Coulomb glass. In such cases, 
the diagrammatic expansion cannot be truncated to any finite order in the interaction, in order to avoid well-known 
divergences associated with the screening processestJ. The simplest consistent treatment has to sum-up all the "chain" 
diagrams, leading to the Debye-Huckel approximation^], which is also known as the random-phase approximation in 
the electronic contexto. This class of diagrams is, in fact, equivalent to the class of "loop" diagrams in our expansion 
of the free energy, which provides the leading contributions in the limit of large coordination. Thus, to address the 
behavior of the Coulomb glass model, one should extend our calculation to sum up all the loop diagrams, which can 
straightforwardly be done even within our formulation. The resulting theory should be capable of addressing the 
interplay of screening effects and glassy freezing, a topic of great relevance for disordered electronic systems. 

D. Weak random fields and self-generated glass 

Another instance where terms not contained in the simplest J 2 theory may be important is the limit of weak 
randomness. On general grounds, one there expects the system at low temperature to assume some uniform order. 
In many cases the corresponding transition has a first-order character, so the upon rapid cooling the system may 
remain trapped in a meta-stable state and undergo glassy freezing. This process is believed to be even more likely 
in presence of competing uniform interactions, which typically can depress the uniform ordering—down to very low 
temperatures. A model for this behavior has been proposed a long time agca, based on earlier workE3 that emphasized 
the relation between mode-coupling theories^ of supercooled liquids and a special class of infinite range spin-glass 
models displaying a first order glass transition scenario. The possibility of glassy freezing in absence of randomness has 
recently attracted renewed attention, and several studiesalU have concentrated on uniformly frustrated infinite-range 
models where these processes can be studied in detail. However, a more general approach would be even more useful, 
where one could examine the competition between uniform and glassy orderings in a controlled scheme, and which 
could be applied to models with realistic interactions and lattice geometries. 

In principle, these questions can be addressed within our approach by examining the fate of the glassy phase in 
the limit of weak random fields, an idea that was introduced a long time agcH In the simplest J 2 theory, and 
outside the ferromagnetic phase our model maps to the SK model in presence of random fields, leading to the glass 
transition line given in Fig. 1(b). As the random fields are reduced, the glassy phase is enhanced, but for sufficiently 
weak randomness, glassy freezing is preempted by uniform ferromagnetic ordering. On the other hand, if we restrict 
our attention to the non-magnetic (m = 0) solution, then the glass transition line can be extended to Hrf = 0, 
leading to Tg{Hrf = 0) = J j \fz. In this way, our formulation may be considered the simplest approach that 
can lead to glassy behavior in absence of randomness. However, the prediction of this lowest order approximation 
cannot be considered as reliable for weak random fields, since corrections to any order in J make contributions of 
the comparable magnitude, even in the limit of large coordination. In fact, using expressions that we have obtained 
within the J 2 theory, it is not difficult to compute the resulting correction to the glass transition temperature, which 
remains of order 1/y/z, but with an increased prefactor. The details will not be elaborated, since stopping at any 
finite order in J is clearly not sufficient. The enhancement of the glass phase in the limit of weak random fields due 
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to these fluctuation corrections may indicate the possibility that, once all the leading corrections are retained, the 
glass transition would precede any uniform ordering as suggested by the SCS theory of Mezard and YoungO. If this 
is correct, it would indicate that the convergence of the 1/z expansion is not uniform as a function of the random 
field strength, since Tg{Hb,f — > oo) ~ 1/y/z, but Tq^Hrf — * 0) ~ 0(1). In that case, the correct mean-field theory 
should not be formulated by performing a z-dependent rescaling of the interactions and then letting z — ► oo. Instead, 
the formulation should retain all the leading 1/z corrections, in finite dimensions. It is interesting to note that recent 
work of Lopatin and Ioffa3 is closely related to the approach that we propose, since it singles out all the leading 1/d 
corrections for a specific uniformly frustrated model, in the limit of large dimensions. In this particular model, the 
interactions are sufficiently frustrated, precluding any uniform ordering, and resulting in a particularly simple large 
coordination limit. In more general cases, competition between uniform and glassy ordering and the possibilities of 
having either first or second order transitions should be considered and may be determined by the details of the 
interactions or the presence of disorder. 

VI. CONCLUSIONS 

In this paper, we have presented a systematic approach that can incorporate short-range fluctuation corrections 
to the standard Bragg- Williams theory of the random-field Ising model. We have shown that if the correct order 
parameters are introduced, corrections to even the lowest nontrivial order immediately result in the appearance 
of a glassy phase for sufficiently strong randomness. This low-order treatment is shown to be sufficient for large 
randomness, where it provides the leading corrections in the limit of large coordination. The structure of the resulting 
glassy phase is very similar to that found in familiar infinite range spin glass models, and is characterized by universal 
behavior emerging from the self-organized criticality of the ground state. 

The major puzzle that remains to be resolved is the extent to which the application of these mean-field ideas is 
relevant to the low temperature behavior of low dimensional systems with short range interactions. An alternative 
approach, based on droplet arguments^ presents a very different scenario, particularly in situations where external 
fields, either uniform or random, explicitly breaks the symmetry of the Hamiltonian. In this instance, droplet argu- 
ments would preclude the existence otany finite-temperature glass transition, in contrast to the mean-field predictions. 
In addition, recent numerical results^ on d = 3 RFIM have also been been used to argue against the existence of a 
finite-temperature transition in a field. In this context, it is worth noting that self-organizedxriticality is not found in 
recent studiescl of hysteresis and avalanche behavior of RFIM with short-range interactions!] in low dimensions. For 
these models, although hysteresis behavior is present, the distribution of avalanche sizes is bounded, and criticality 
is found only by fine-tuning the parameters of the system to a particular point of the phase diagram. Similar results 
have been obtained in studies that have examined the sensitivity of the d = 1, 2 RFIM to small random perturbations 
of the quenched disordered. 

In our opinion, more general emergence of self-organized criticality similar to that found in mean-field glassy models 
most likely requires-the existence of longer-range spin-spin interactions and/or high spatial dimensions. On the other 
hand, experimentscj measuring the Barkhausen noise on several "hard magnets" have indicated that power-law 
distribution of avalanche sizes and avalanche times, consistent with self-organized criticality. Such behavior may be a 
result of the fact that in such systems the dominant interactions have a dipolar, thus longer-range character, bringing 
the behavior of these materials closer to the predictions of mean-field glassy models. These features may also be of 
particular importance in applications of the RFIM to disordered electronic systems and the related physics of the 
Coulomb glass behavior. 

In this work, have also outlined how our theory could be extended to examine models with either longer-range 
interactions or the limit of weak random fields, which is of particular importance to the long-puzzling_queatipn of 
the self-generated glassiness in uniform systems. Our theory is closely related to other recent approachesEj'oEH3 that 
address the emergence of glassy phases on a mean-field level. These theories taken as a whole appear to shed light on 
a number of experimentally relevant systems, and present a fairly complete and consistent picture of glassy behavior. 
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APPENDIX A: CANCELLATIONS DUE TO UP-DOWN SYMMETRY 

In investigating the glassy phase of the RFIM outside the FM phase, and in zero uniform field, we can make use of 
the fact that after averaging, the system respects up-down symmetry. As a result, a number of terms vanish, so that 
the expressions simplify. In the following, we discuss these cancellations in some detail. 

1. Moments of spins on different sites 

The term (</>) contains expectation values of a product of spins on different sites. Since the averages are computed 
at s = 0, such terms very generally factor out, so that we can write 

(St 1 • ■ • S? k ») = (Sg)o ■ ■ ■ (S? k k )o = m fc , (Al) 

where m is the magnetization. Outside the FM phase, and at H — 0, the magnetization and thus all such moments 
vanish. 

2. Evaluation of ((f)) and (\)o 

An immediate consequence of the above factorization property is the fact that (<p) — 0, since this expression 
contains products of two spins on different sites. The expression for \ contains a derivative of the source field of the 
form -§-£,"■ To compute this derivative, we use the fact that due to the definition of the Legendre transform, we can 



3£M 

write 



& = tt* (/3G) (A2) 



To calculate (x)o, we need to compute the derivative of £? (e), evaluated at e = 0, and we get 



dt; 



a0 



Be 

We therefore conclude that (x)o = as well 



l) m 



£=o dq fde 



Wo 



e=0 ^f 



(A3) 



3. Terms with odd powers of <f> 

In evaluating higher order terms in the e-expansion, terms of the form (4> p tp(x' \x"))o appear, where i/j(x'iX") i s 
an arbitrary polynomial function of \' and %", and p is an odd number. To evaluate such terms, we note that as 
before, spin moments on different sites factor out, but we still have to compute nontrivial spin moments of the form 
(S" 11 ■ ■ ■ S°' r ) , with r having the same parity as p, i .e r is odd. Here, we have used the fact that quantities x' = -§^X 

and x" = TpzX are quadratic in local variables, i. e. contain products of the form SfSf. In addition, the considered 
moments will have an odd number of spins only if the considered lattices have no odd-membered rings, such as found 
for example on a triangular lattice. If any of the replica indices coincide, then (S") 2 = 1, and an even number of 
spins drop out, but the remaining expression still has the form 

M S = (S^---S^) , ai <a 2 <---a s . (A4) 

Using well-know properties of replicasc3, it is readily seen that 



M a = (Si)' = m s = 0. (A5) 

Therefore the expression of the considered form vanish as well for lattices with no odd-membered rings. 
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APPENDIX B: VARIATIONS OF THE GIBBS FREE ENERGY 

A general procedure needed to obtain the equation of state and the RSB stability criterion involves computing the 
variations of the Gibbs free energy with respect to the order parameter q a/3 . In the following we outline how these 
variations can be computed by concentrating only on the leading-order contributions from the "loop" diagram. 

1. Calculation of d a pgint 

The calculation of d a pgi n t to order J 2 is already computed in section III A. here, we compute the J 4 correction. 
Defining the matrix A = I + Q, we can write 34 = -4-TrLA 4 ], and find 

a n -V" d9A dA "< s - 2i (P\ rm\ 

9^34-X, ^ - ~M U. (Bl) 

70 

In the RS limit q a/3 = q, and taking n — > this reduces to the expression of Eq. (41). 

2. Calculation of g" nt 



From Eq. (31), the J 2 contribution reads 

2 
z 



d lf3,jS 92 = -SaPrfS- (B2) 



The J 4 contribution can be calculated using the same procedure as for the first variation, and we find 

48 



d a /3,~{S 94 



Sari 2_^ Ap^A^s + A ai Aps + dps 2_^ A atl A^ 



(B3) 



In the replica symmetric theory, the resulting matrix elements of g" nt are given by 



-2 1 /„2 



2 



P = (SL)afi,afi = £ ~ + \ (j) (2(1 - q 2 ) + 1), 

1 fe 2 \ 2 
Q = {9int) a p, a s = 5 [ — ) ( 2 9(! - l) + <l), 

R=(9intU, 1 8 = U E A Q 2 - (B4) 

Given these matrix elements, one can immediately evaluate the relevant eigenvalueEa of the matrix g'/ nt , which takes 
the form 

Af' = P-2Q + iU^ + |(^) (1-g) 2 . (B5) 

In this expression, we note that similarly as in the computation of the J 2 correction to the RS equation of state, the 
J 4 correction is proportional to 6q = 1 — q, and is therefore down by a factor (J / Hrf) 2 compared to the leading term. 
We again conclude that in the limit of large random fields, to leading order it suffices to retain the J 2 contribution. 

3. Calculation of /" 

The functional / [£] is the free energy of free spins in presence of fields £ Q/3 , and we are interested in computing 
its second variation at the saddle point where S^ 13 = d a 09int- Therefore, this evaluation is almost identical as for 
the SK modeEJ. In addition, since this quantity is evaluated in the RS theory, dapgint assumes the form that we 
have already discussed when we examined the RS equation of state. We conclude that for this quantity as well, to 
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leading order in the limit of large random fields, the argument of /" can be replaced by its J 2 approximation, and we 
straightforwardly obtain Eq. (11). 
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